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ABSTRACT 

Q . The calculation of the transmission power spectrum of QSO's Lya 

, absorption requires two parameters for the normalization: the continuum F^ 
, . and mean transmission e""^. Traditionally, the continuum is obtained by a 

1 ^ . polynomial fitting truncating it at a lower order, and the mean transmission 
Q . is calculated over the entire wavelength range considered. The flux F is then 
^~~' . normalized by Fce~'^ . However, the fluctuations in the transmitted flux are 
y—^ . significantly correlated with the local background flux on scales for which the 
^ , field is intermittent. As a consequence, the normalization of the entire power 
OO . spectrum by an over-all mean transmission e""^ will overlook the effect of the 
r^^ ! fluctuation-background correlation upon the powers. In this paper, we develop 
O , a self-normalization algorithm of the transmission power spectrum based on a 
Q ! multiresolution analysis. This self-normalized power spectrum estimator needs 

r-| ! neither a continuum fitting, nor pre-determining the mean transmission. With 

^! simulated samples, we show that the self-normalization algorithm can perfectly 

2 '. recover the transmission power spectrum from the flux regardless of how the 

c/3 ! continuum varies with wavelength. We also show that the self-normalized power 

. . \ spectrum is also properly normalized by the mean transmission. Moreover, this 

• y-{ ! power spectrum estimator is sensitive to the non-linear behavior of the field. 

r> \ That is, the self-normalized power spectrum estimator can distinguish between 

5^ \ fields with or without the fluctuation-background correlation. This cannot be 

accomplished by the power spectrum with the normalization by an overall 
mean transmission. Applying this analysis to a real data set of ql 700-1-642 
Lya forest, we demonstrate that the proposed power spectrum estimator can 
perform correct normalization, and effectively reveal the correlation between the 
fluctuations and background of the transmitted flux on small scales. Therefore, 
the self-normalized power spectrum would be useful for the discrimination 
among models without the uncertainties caused by free (or fitting) parameters. 
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1. Introduction 

Lya absorption, shortward of Lya emission in QSO spectra, indicates the presence 
of intervening absorbers with neutral hydrogen column densities ranging from about 10^^ 
to 10^^ cm~^. The absorbers with low column densities, e.g. from 10^'^ to 10^^ cm~^, are 
usually called Lya forest. It is generally thought that the low column density absorbers are 
some kind of weakly clustered clouds consisting of photoionized intergalactic gas {e.g. Wolfe 
1991; Bajtlik 1992). This suggests that Lya forests are caused by diffusely distributed 
IGM in pre-collapsed areas of the cosmic mass field (Bi, 1993; Fang et al. 1993; Bi, Ge 
& Fang 1995; Bi & Davidsen 1997.) Observations of the size and velocity dispersion of 
the Lya clouds at high redshift also show that the absorption probably is not caused by 
confined objects at high redshifts (Bechtold et al. 1994; Dinshaw et al. 1994; Fang et al 
1996; Crotts & Fang 1998.) With this picture, the baryonic matter distribution is almost 
point-by-point proportional to the dark matter distribution on all scales larger than the 
IGM's Jeans length, i.e. the Lya forests would be good tracers of the underlying dark 
matter distribution. 

Thus, the power spectrum of QSO Lya transmitted flux can be used to estimate the 
power spectrum of the underlying mass fleld, and then be used to constrain cosmological 
parameters (Croft et al 1999; McDonald et al 1999; Hui 1999; Feng & Fang 2000.) A 
key step in this approach is to compare the power spectrum of observed transmitted flux 
fluctuations with model-predicted power spectrum. One uncertainty in the power spectrum 
determination of the real data is from the normalization of the power spectrum. Therefore, 
in order to have an effective confrontation between the observed and theoretical power 
spectrum of Lya forests, it is necessary to develop a proper algorithm for the normalization 
of the power spectrum. This is the goal of this paper. 

The observed flux of a QSO absorption spectrum is given by F{X) = Fc(A)e~^, where 
-Fc(A) is the continuum, e^^*^'^-' the transmission, and r the optical depth. The normalized 
power spectrum of transmission is the power spectrum of the transmission flux fluctuations 
5(A), deflned as 

^ ' F(A) ^ ^ 



That is, the transmission power spectrum is normalized by the mean flux F{X) = Fc{X)e'^^'^\ 
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In other words, the normahzation of the transmission power spectrum is determined by two 
factors: the continuum Fc{X) and the mean transmission e^('*'). 

Traditionally, the continuum is needed to be determined before the power spectrum 
calculation. Usually the continuum is obtained by a fitting of polynomial or its variants. 
Assuming that the continuum fiuctuates slowly, the polynomial or its variants are truncated 
at relatively low orders (e.g. Croft et al 2000; Hui et al 2000). The pre-assumed polynomial 
or other function, and the subsequent truncation may lead to uncertainty of the power 
spectrum. 

Another source of uncertainty of the transmission power spectrum is the mean 
transmission e~^ normalization. The mean transmission is calculated by averaging the flux 
over the entire wavelength range considered, and the power spectrum is normalized by 
this mean transmission for all scales. This implicitly assumes that there is no correlation 
between the transmitted flux fluctuations and the mean flux. This assumption is true for a 
gaussian fleld, but may not be so for a non-linearly evolved fleld. 

In fact, the fluctuations at position A are correlated with the background at the same 
position A. Recent flndings that the transmitted flux of Lya forests exhibits intermittent 
behavior (Jamkhedkar, Zhan &; Fang 2000) clarifles this point. That is, the transmitted 
flux shows prominent spiky feature fluctuations on small scales. The transmitted flux 
consists of rare but strong density fluctuations randomly scattered in space with very low 
fluctuations in between. In this case, the power of the transmission fluctuations is mainly 
dominated by the spikes. On the other hand, the transmission e""^ is low at the spikes. That 
is, the transmission fluctuations are anti-correlated with transmission. As a consequence, 
the power would be underestimated if the power spectrum is normalized by the mean 
transmission e~^ over the entire wavelength range. Since the spiky features are stronger on 
smaller scales, the normalization by an over-all mean transmission e~'^ or by a fllling factor 
with a scale-independent flux threshold (Croft et al 2000) will cause an underestimation of 
power on small scales. 

Recently, we have developed a power spectrum estimator with a multiresolution analysis 
based on the Discrete Wavelet Transform (DWT). The DWT power spectrum estimator 
is found to be very useful for the recovery of the initially linear power spectrum (Feng & 
Fang 2000; Pando, Feng & Fang, 2001). In this paper, we show that the DWT algorithm is 
also very useful to detect the power spectrum of non-linear fleld, like an intermittent fleld. 
The DWT algorithm can effectively reduce the above-mentioned uncertainties due to free 
parameters used for normalization. We will show that the normalization of a DWT power 
spectrum does not rely on a continuum fltting and the mean transmission. Moreover, the 
power spectrum given by this estimator is sensitive to the correlation between the flux 
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fluctuations and the background flux. That is, the power spectrum can be employed to 
distinguish among the fields with and without intermittency. Therefore, it would be useful 
for discrimination among models of the Lya forests. 

The paper will be organized as follows. §2 introduces briefiy the Discrete Wavelet 
Transform (DWT) analysis of the fiux of QSO absorption spectrum. §3 presents the 
self-normalization algorithm of the transmission power spectrum. It doesn't need either 
a continuum fitting, or a calculation of the mean transmission. In §4, we test the self 
normalization algorithm. We show that the self-normalization algorithm can effectively 
perform the normalization due to either continuum or mean transmission. §5 will 
demonstrate that the DWT self-normalized power spectrum is useful to detect the 
intermittent behavior of the field. Finally, the conclusions and discussions will be presented 
in 56. 



2. The DWT analysis of the QSO Lya forest flux 

2.1. Need for a space-scale decomposition 

We rewrite eq.(l) as 

F{X)=F{X} + F{X}S{X). (2) 

Our purpose is to estimate the power spectrum of the fiux fiuctuations 6{X) from the 
observed fiux F{X). Therefore, eq.(l) requires one to decompose the observed fiux F{X) 
into two terms: the first one is the background F{X), which does not contain information 
of the fiuctuations considered, while the second term contains all this information. If 
the background F{X) is A dependent and correlated with the fiuctuations S{X), the 
decomposition eq.(l) apparently cannot be done by truncating at a priori "relative low 
orders" . 

We will solve this problem by a scale-by- scale analysis, without introducing new 
parameters. In terms of a scale-by-scale analysis, eq.(l) means that to detect the power of 
the fiux fiuctuations on the scale r, all components of F{X) on scales larger than r play the 
role of a background. Therefore, to determine power on the scale r, one can decompose the 
observed F{X) into two terms: the first doesn't contain any information on scales larger 
than r, while the second contains all information on scales equal to or less than r. 

This decomposition refers to both the position A and the scale r, and therefore, we 
need a scale-space decomposition of the transmitted fiux. 

On the other hand, the calculation of power spectrum essentially is a decomposition 
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of the flux into scale domains. Therefore, it is possible to do the decomposition of eq.(2) 
by the same scale-space decomposition as that used for measuring the power spectrum. In 
other words, the estimation of the normalization background and the calculation of power 
spectrum can be accomplished simultaneously. Once the orthonormal bases for power 
spectrum estimation are given, the term F{X) is uniquely determined without free (or 
fitting) parameters. 

The Fourier power spectrum is not convenient for this purpose, as the bases of 
Fourier transform are not localized in physical space, and they don't yield a scale-space 
decomposition. We will use the DWT, whose bases are localized both in scale and position 
space (e.g. Mallat, 1989a,b,c; Meyer, 1992; Daubechies, 1992; and references therein, and 
for physical applications, refer to Fang & Thews, 1998.) 



2.2. Expansion by scaling functions 

A sample of QSO Lya absorption spectrum gives a list of the flux F{\) observed at 
discrete wavelength Aj, i = l,2...n. Since Aj corresponds to spatial position Xi, or redshift 
Zi, one can define a spatial distribution of the fiux by 

n— 1 1 

Fi^) = E oinK) + F{K^,)][e{x - A,) - 9{x - A,+i)], (3) 

where the step function 6{x) is equal to 1, for a; > 0, and for x < 0. The spatial range L 
corresponds to the wavelengths from Xmin to Xmax- 

In the DWT analysis, the space L is chopped into 2^ segments labelled by 
/ = 0, 1, ...2-' — 1. Each of these segments has a size of L/2^ . The index j is an integer. 
Therefore, the index j stands for scale L/2^ , while / is for position, or the spatial range 

lL/2^ <x <{l + l)L/T. 

We first introduce the scaling functions for the Haar wavelets. They are the top-hat 
window functions defined by 




(4) 
otherwise. 



where the superscript H stands for Haar, and the factor J2^ / L ensures normalization, i.e. 

f 

Jo 



<pf,i{x)<pf^Ax)dx = 6^„ (5) 
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where 5^ is Kronecker delta function. 



The scaling function, (f)fi{x) is actually a window on a resolution scale L/2^ at the 
position L12~^ < a; < L[l + 1)2^-'. But it is not normalized like a window function W{x) 
which satisfies 

''w{x)dx = l. (6) 



/ 



Nevertheless the mean flux in the spatial range L12 ^ < x < L{1 + 1)2 ■'is proportional to 

ef,= ("^ F{x)c^Ux)dx. (7) 



The number e^i is called the scaling function coefficient (SFC). 
Using SFCs, one can construct the flux as 



F'^^) = E 4i<pU^)- 



1=0 



F^{x) is the flux F{x) smoothed on scale L/2^ (or simply the scale j). A higher value for 
j corresponds to smaller scales and vice versa. For a given sample with resolution 6X, the 
original flux can be expressed as F{x) = F'^ , where J is given by the integer of number 



2.3. Expansion by wavelets 

F^{x) contains less information than F{x), because information on scales > j + 1 
(i.e. smaller scales) has been smoothed out. It would be nice not to lose any information 
during the smoothing process. This can be accomplished if the differences, F^{x) — F^'^{x), 
between the smoothed distributions on successive scales are retained. In other words, only 
if we are able to retain all these differences, this scheme will not lose any information. 

To calculate the differences, we define the difference function, or wavelet, as 

.2,. 1/2 [ 1 for0<^<l/2 
tlj^(r]) = y < -1 for 1/2 < 7] < 1 (9) 

I otherwise. 

This is the basic Haar wavelet. One can then construct a set of wavelets "ipuix) by dilating 
and translating eq. (9) as 

^5(x)=^^(2V^-/). (10) 
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The Haar wavelets are orthonormal with respect to both indices j and /, i.e. 

For a given j, ip^i{x) is also orthogonal to the scahng functions 4'y,i{.x) for f < j, i.e. 

(t>lv{x)i,fi{x)dx = Q, for j'<j. (12) 

-^ J. 



(13) 



From eqs.(4) and (9), we have 

0S/(a^) =0f-i,/(a;)+^f_i/x), 

Thus, the difference F^{x) — F^~^{x) is given by 

F^{x) - F^'\x) = ' E '^Ti,/<-u(^)' (14) 

1=0 

where ej^_^ ; are called the wavelet function coefficients (WFC) given by 

~ef,= j F{x)Mx)dx. (15) 



Using the relation (14) repeatedly, we have 

F^{x) = F\x) + x: '^Y:~4,i^f'A^)- (16) 

j,'=0 1=0 

This is an expansion of the flux F^{x) with respect to the basis il)fi{x), and F^{x) is the 
mean of F{x) on the entire range L. Therefore, the flux F{x) can be can be expressed as 

F{x) = F\x) (17) 

j'=j 1=0 
j'=0 1=0 

The Haar wavelet provides a clear picture of the DWT decomposition, and it is also 
easy for numerical work. However, the Haar wavelet is discontinuous in real space, and 



therefore, it is not well behaved in scale space. For our work, the most important properties 
of the basis for the scale-space decomposition are 1.) orthogonality, 2.) completeness, 
and 3.) locality in both scale and physical spaces. Therefore, all wavelets with compactly 
supported basis will produce similar results. Among the compactly supported orthogonal 
bases, the Daubechies 4 (D4) is easy for numerical calculation. We will use wavelet D4. 

For the D4, the basic orthonormal eqs.(5), (11) and (12) still hold after replacing 0^ 
and ip^i with D4 scaling function (pji and D4 wavelet ipji. The expansion eqs.(17) are also 
valid using (pji and ipji. 



2.4. The DWT power spectrum of flux without normalization 

Form eq.(17), the Parseval theorem of the DWT gives 



J[F{x) - F^dx = J:T. \e 



~F|2 
j=0 1=0 



Therefore, the power of the mode (j, /) is |ef;p. Thus, the power spectrum of the flux F{x) 
is given by (Pando & Fang 1998; Fang & Feng, 2000; Yang et al. 2001) 

pf = {r4/) (19) 

The ensemble average of the power (jejl^jp) should be /-independent for a homogeneous field. 
If the "fair sample hypothesis" is true (Peebles 1980), the ensemble average can be replaced 
by a spatial average, we then have 

Pf = ^E\~4/- (20) 

This is the DWT power spectrum of the flux F{X). It is a band-averaged Fourier power 
spectrum as 

-| oo 

Pj = ^ E mnmmn), (21) 

where P{n) is the Fourier power spectrum with wavenumber k = 27Tn/L, and ■?/'(n) is the 
Fourier transform of the basic wavelet iIj{x). This relation has been conflrmed with 1-D 
(Lya) and 3-D (N-body) numerical samples (Feng & Fang 2000; Yang et al. 2001.) 

The DWT algorithm eq.(21) have been successfully employed for the power spectrum 
reconstruction (Feng & Fang 2000; Pando, Feng & Fang 2001). It can recover the initial 
linear power spectrum on small scales as well as large scale. For the reconstruction, the 
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normalization is simple, as the field is gaussian. For non-gaussian fields, the normalization 
will no longer be trivial, as the non-gaussianity will affect the normalization differently for 
different algorithm. 



3. Algorithm for the normalized DWT power spectrum 
3.1. The DWT power spectrum from a Poisson sampling 

Using photon numbers N{x) oc F{x), eq.(2) can be rewritten as 



N{x) = N{x) [1 + 6{x)] + Ne{x). (22) 



where N{x) = Nc{x)e'~'^^^^ is the mean photon number at x (wavelength or redshift). The 
term Ne{x) is given by noise, which includes the background sky, dark current and the 
instrumental readout noise (see Appendix A). 

By definition, the mean photon number (or mean fiux) is the background of the 
fluctuation 5{x). It corresponds to the selection function in the problem of galaxy 
distribution. Because N[x) or F{x) are not constant in the entire spatial range, the 
background can only be defined scale- by-scale. In other words, to measure the power of 5{x) 
on a given scale j, the background N{x) and F{x) are given by the mean photon number 
and mean fiux at x when all fluctuations on scales smaller than the given scale (> j) are 
absent, i.e. 



N{x) = N{l).[l + d{x),] + Ae(x). (23) 

and 



F{x) = F{l)^[l + 6{x),]+e{x). (24) 

where 5{x)j contains all information about the transmission fluctuations on scales > j. 
N{x)- and F{1)- are the mean photon number and mean flux in the spatial range 
lL/2^ < X < {I + l)L/2\ respectively, i.e. the local mean photon number and mean 
flux when the mass field clustering on scales > j is absent. The error term in eq.(24) is 

e(x) = Ne{x)F{x)/N{x) as N{x) oc F{x). 

Due to the discreteness of photons, the observed photon number N{x) should be 
considered as a sampling of the random field of eqs.(22) or (23). Our purpose is to estimate 
the power spectrum of the random field S{x), which describes the fluctuations of the 
transmission. If the sampling is Poissonian, eq.(22) yields 
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-5''ix-x') 



'N,{x)N,{x') 



,N{x)l \N{x)N{x')l 

The derivation of eq.(25) can be found in Appendix B which describes both Poisson and 
non-Poisson samphng. Considering F{x) oc N{x), the second and fourth terms on the r.h.s 
eq.(25) can be rewritten as 

F{x)F{x') \ 



{6{x)6{x')) = -1 + 
-6''(x-x') 



\F{x) F{x')/ 

1 \ / e{x)e{x') 



(26) 



N{x) 



Fix) Fix' 



Projecting eq.(26) onto the DWT bases ipji{x)ipji{x'), the l.h.s. gives the DWT power 
spectrum of S{x), i.e. the normalized DWT power spectrum of the flux fluctuations, 



^ = (4)- 



(27) 



The first term on the r.h.s. of eq.(26) disappears, because 4'ji{x) is admissible, i.e. 
J ipji{x)dx = 0. Thus, eq.(26) yields 

2, 



Pi 






I 1 



\n{i), 



~E 



(28) 



where e^ = j' ei^x)ipj^i{x)dx. 

Using again the "fair sample hypothesis" , we have 



P^ 



2^ 



2J-1 

E 

1=0 
2J-1 






(29) 






2^ 



a N{1)^ 

In calculating the error term, we used 

(e(A)) = 0, 



2J-1 

E 

1=0 



2^ f^ [F(/)P 



a'^{x)ilj^i{x)dx. 



(e(A)e(A')) = ^^(A)^^,^ 



(30) 



The standard derivation o"(A) can be found from the data set given by observers. It is also 
reasonable to assume that the noise is uncorrelated with 5{x). 

Eq.(29) is the basic estimator of the normalized DWT power spectrum of Lya 
transmitted flux. The first term on the r.h.s. of eq.(29) is the power spectrum of the 
transmission. The second and third terms are the corrections due to shot noise and the 
noise e(A), respectively. 
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3.2. The estimation of the mean flux 



Now we turn to the problem of estimating the mean flux F{1)-. The deflnition of F{1)- 
given by eq.(24) requires a decomposition of F{x), in which the first term contains only 
information on scales larger than L/2^ , but nothing on scales less than or equal to L/2^ , 
and all the information of flux fluctuations on scales equal to or less than L/2^ is contained 
in the second term. 

This decomposition is already given by second line of eq.(17), in which F^{x) is the 
flux at X if all the fluctuations on scales less than or equal to L/2^ are absent. Thus, to 
measure the power on scale j, one can identify the mean flux F{1)- with F^{x). Obviously, 
F^{x) is not constant, except within each small segment lL/2^ < x < {I + l)L/2\ as all 
fluctuations on scales > j are smoothed out. Therefore, we have 

m,=F'i^)='i:4M^), (31) 

1=0 

where, for a given /, x is given by lL/2^ < x < {I + 1)L/2K 



Eq.(31) shows that the mean flux F(l)- at position / can be represented by the SFC ej; 
for scale j. However, as mentioned in §2.1, the SFC ej^ is proportional to, but not equal 
to the mean flux at position /. To flnd the proportionality constant, we use the so-called 
"partition of unity" of wavelets given (Daubechies 1992) 

oo 

^0(r^-/) = l, (32) 

Z=— oo 

or 

, T .1/2 2^-1 

1^) E <^j^(^) = 1' (33) 

where (t)fi{x) is the so-called periodized scaling function (Fang & Feng 2000). With 0j^(a;), 
eq.(7) can be rewritten as 

eli= 1"^ F{x)cPl,{x)dx. (34) 



Substituting eq.(3) into eq.(34), we have 

"-1 1 r\ 



^U 



Y: -[F{X.) + F(A.+i)] / '^' <j>f^,{x)dx. (35) 

Thus, with the "partition of unity" eq.(32), eq.(35) yields 

2-'-l . r s 1/2 n-l -r 

E Uy <' = ^ 9f^(^^) + ^(^^+i)](Am - A.). (36) 

1=0 ^^^ i=l ^ 
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This result is clear. The r.li.s. of eq.(36) is the integral flux of the entire spectrum, while the 
l.h.s. is the integral flux of the smoothed spectrum F^{x). That is, the smoothed spectrum 
F^{x) always has the same integral flux as the original spectrum. F^ is only a reassignment 
of the original flux F{x) from the distribution F{Xi) over grids i = 1.. n to new distribution 
^jL/2^efi over grids l = 0..2^ -1. 

The quantity (L/2^y/2gF jg ^j^g g^^ j^ ^j^g spatial range lL/2^ to (/ + l)L/2^. Thus, 
the mean flux in this spatial range is 



2^ 



1/2 



F(l). = {j] 4- (37) 



Although the mean flux F(l)- is different on different scales, the average of F{1)- over / is 
independent of j, i.e. 

This allows comparison of powers on various scales. In other words, the power, Pj, on 
different scales j is properly normalized. 



From the definition of the mean flux, F{1)- can be written as 



F{l)^=F^{lUe-r%. (39) 

where Fc{l)j is the continuum at the cell (j, /). 

If the continuum can be approximated as a constant, i.e. Fc{l)j = Fc, eq.(38) gives 

^'ll^^j = ^yE ^[^(A.) +F(A.+i)](A,_,i - A,). (40) 



Eq.(40) means that the average of the mean transmission e~^(')j over positions / is 
independent of j, i.e. it is a constant. However, even in this case, it doesn't mean that the 
power spectrum can be normalized by a constant mean transmission (§4). 



3.3. The estimator of normalized power spectrum of transmission 

Substituting the mean flux eq.(37) into eq.(29), we have the transmission power 
spectrum estimator as 



1=0 



4 



T 2J-lr2F-|^ -I 2J-1 1 r 2^-1 2(1\ 



2. g Nil), 22. g [ej^]2 ' (41) 
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where 



a^{l) 



a 



x)ijj'^i{x)dx. 



(42) 



It is the mean of a"^(/) at the position /. 



In the algorithm of eq.(41), the fluctuation amplitudes e^i (WFCs) at position / are 
normalized by the mean flux ej^ (SFCs) at the same position /. That is, the normalization 
coefficients (SFCs) and the fluctuation amplitudes (WFCs) are obtained from the same 
DWT decomposition eq.(17). We call it the self-normalized power spectrum estimator. 

Let consider the case of a constant Fc. Eq.(29) and (39) give 



Pi 



C. 



Fce--«, 



.^(0., 



~E 



F.e-(% 



(43) 



A constant F^. actually doesn't affect the power Pj, because the power is given by the ratio 
efi/efi [eq.(40)], and both the ifi and e^ are proportional to Fc [eqs.(7) and (34)]. 



If there is no correlation between e^i and e '^(^)j, eq.(43) gives 



Pj 



|e--(0^f, 



4i 



I 1 



\N{1), 



~E 



n2. 



i^ce--«,. 



(44) 



From eq.(40), the factor (l/|e~'^(')jp) is approximately independent of j. That is, the 
normalization of power spectrum (44) is given by a constant mean transmission. If the 
WFCs eji are correlated with e^^(')j or the SFCs ej';, the power spectrum cannot be 
normalized by a constant mean transmission e~'^, even when the continuum is constant. 



4. Test of the self-normalization algorithm 

In this section, we show that the power spectrum given by eq.(41) is properly 
normalized. It doesn't need a continuum-fitting, or a pre-calculation of the mean 
transmission. 



4.1. Simulation samples of the lognormal model 



Our purpose is only to demonstrate the normalization, but not model discrimination or 
cosmological parameter determination. We can use the semi-analytical model, the so-called 
lognormal model developed by Bi & Davidsen (1997). The lognormal distribution is useful 
for testing algorithm on intermittent fields, as lognormal field is typically intermittent. 
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Briefly, this simulation consists of three steps. First, we generate a reahzation of the 
primordial (linear) baryonic mass distribution of size 200 h~^ Mpc in comoving space, 
centered at a typical redshift, say z = 2.4 in fi = 1 cold dark matter universe. The COBE 
data is used to normalize the initial spectrum which has the soften parameter F = 0.3. The 
Hubble constant is taken to be 50 km s~^ Mpc~^. In the Fourier space, the distribution 
of baryonic matter differs from that of dark matter due to smoothing on the Jeans length. 
Second, the non-linear baryonic density is calculated using the lognormal transformation 
from the linear density; and the neutral baryonic density is calculated assuming an UV 
ionization flux of Jq = 1.0 x 10"^^ erg cm^^ s^^ Hz^^ sr^^. The temperature of the gas 
is assumed to have a mean of 2.0 x 10^ K following an adiabatic equation of state with 
the polytropic index 4/3. Because of the UV radiation, a minimal temperature should 
be introduced. We take it to be 1.12 x 10^ K. Lastly, the absorption of the Lya photons 
by hydrogen is convolved by a Voigt profile, and the whole optical depth is calculated by 
summing over all the pixels. Each pixel has the size of 0.0156A and the total number of 
pixels is 32768. 

To fit with the observed data of medium resolution, we convolve the theoretical 
spectrum by an instrumental point-spread-function (PSF) of the typical resolution 41 km 
s~^. Alternatively, one can simulate a medium resolution spectrum by using a higher gas 
temperature but without the instrument PSF. We will use the temperature 1.0 x 10^ K to 
mimic the instrumental effect. With this simulation, we have the transmission e^'^^'^\ which 
is shown in the top panel of Fig. 1. 

We then assume the CCD counts Nc = 2000 per wavelength if the fiuctuations are 
absent. Three continua Fc are considered: 1) flat continuum, Fc = N^ 2) Fc = Nc{\/\o)^, 
where Aq is 4133 A(redshift 2.4), 3) Fc = Nc[l + dsm[2n{\ — \min)/i^max — ^min)], and 
A„i„ = 3890.49 A, X^ax = 4399.53 A. The amplitudes d is taken to be 0.1, 0.2, 0.3, 0.4 
and 0.5. The observed photons are a Poisson sampling of a random fleld with mean Fce~'^. 
For each given Fc, hundred realizations of the Poison sampling are produced. Finally, 
gaussian noise with zero mean counts and standard derivation 50 (photon number), which 
is independent from the Poisson sampling, is added to each pixel. Some results are shown 
in Fig. 1. 



4.2. Continuum-independent power spectrum 

Using the estimator in eq.(41), we calculate the DWT power spectra from the simulated 
flux with various continua. The results are plotted in Fig. 2. It shows that the DWT 
power spectra for different continua are exactly the same on scales j > 3 (or length scale 
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2^^~^ X 2 h^^ kpc). The dispersion of the power spectra over 100 reahzations is very smaU. 
That is, the power spectrum given by estimator (41) is continuum-independent. That it, 
the self-normahzation algorithm can produce the power spectrum correctly normalized by 
the continuum, but without a continuum fitting. 



4.3. Normalization of mean transmission 

To test the mean transmission normalization of estimator eq.(41), we calculate the 
so-called unnormalized DWT power spectrum of the transmission, i.e. the power spectrum 
of continuum normalized flux 

(45) 



-.(A) _ F(X}_ 



Similar to eq.(20), the unnormalized DWT power spectrum is given by 



^' = ^EI^1/ (46) 

^ 1=0 



where the WFCs are calculated by 



'h = Jjr^^M^)dx. (47) 

The result of P* for the lognormal sample is plotted in Fig. 2. 
By definition eq.(l), we have 

p-rW 

HX) = p^ - 1 (48) 

The second term on the r.h.s. of eq.(48) does not contribute to the DWT power spectrum, 
as ipji{x) is admissible. Therefore, the power spectrum P* should be transformed to Pj by 
dividing a normalization factor [e""^]^, or 



e 



-rW^ 



[p;/p,f/i (49) 



One can test the normalization of mean transmission by comparing the value of e^^(^)j given 
by eq.(49)] with the value e~'^ used in the simulation. For the lognormal simulation samples, 
the result is shown in Fig. 3, in which three samples with different e~^ are analyzed. It 
shows that e^'^^^^j on small j (large scales) is exactly the same as the simulation used e~'^. 
That is, the estimator eq.(41) is able to produce the power spectrum, which is already 
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properly normalized by the mean transmission. Thus, the estimator eq.(41) doesn't need 
one to pre-calculate the mean transmission or filhng factor for the normahzation. 



However, on small scales, e~'^^^'>j doesn't equal to e^^. Therefore, on small scale, the 
estimator eq.(41) is not the same as the traditional normalization by a mean transmission. 
This difference will be analyzed in next section. 



5. Normalization of power spectrum of a non-linear field 
5.1. Effective mean transmission 



We call e^'^^^^j the effective mean transmission. Fig. 3 shows that the effective mean 
transmission e~'^^^^ j is scale-dependent. It means that the power spectrum given by eq.(41) 
actually is normalized by a scale-dependent factor e~^('*')j, not by a mean transmission. 
Generally, e~'^('*')j is lower than the mean transmission e~^ on small scales. 

The scale dependence of the effective mean transmission is due to the anti-correlation 
between the background e~'^'^'^) and the fluctuation 5{X). For a lognormal sample, this 
anti-correlation is shown in Fig. 4. The mean value of the WFCs e^ (fluctuation) is 
larger for smaller SFCs e|J (background), and vice versa. On the other hand, the power 
is dominated by large WFCs, corresponding to smaller SFCs. Thus the normalization by 
an over-all averaged SFCs, or mean transmission e~'^ will produce lower power than the 
power spectrum given by eq.(41). This anti-correlation comes from the spiky structures of 
an intermittent field. It is mostly significant on scale j = 8 and 9 (100 - 250 h~^ kpc). On 
scale j = 10 (or 50 h~^ kpc), the anti-correlation is weak due to noise. 

On the other hand, the field on large scales basically is gaussian, i.e. no correlation 
between the WFCs ej^ and the SFCs ej^. In this case eq.(44) shows that the normalization 
can be done by a constant mean transmission e^'^(^). This is why the mean transmission 
can be recovered by the ratio eq.(49) on large scales [Fig. 3]. 



5.2. Test by pseudo-hydro simulation samples 

In the lognormal model, the non-linear feature of the field is induced by the lognormal 
transform. To further test the normalization of effective mean transmission, we use 
simulation samples of the Lya forests produced with the so called "pseudo- hydro" technique 
(Croft et al. 1998). In this approach, the non-linear density and velocity field of underlying 
dark matter were obtained by evolving a particle mesh (PM) simulation for a specified 
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cosmological model. The gas density and temperature were then computed using the simple 
scaling relation inferred from full hydrodynamic simulations (Hui & Gnedin 1997). 

The parameters of the simulation are taken to match the Keck HIRES spectra described 
in Croft et al (2000). First, the PM simulations were performed by evolving 128'^ dark 
matter particles in a periodic box with a 128^ grid. For this study, the cosmological model 
is taken to be the low density fiat model (LCDM), which was specified by the density 
parameter Qq = 0.3, the cosmological constant Q\ = 0.7, the Hubble constant h = 0.7 and 
the baryon density parameter fit = 0.0125/i~^. The physical size of box was determined 
correspondingly by 512 pixels around the median redshift Zmed of the sample, and a total of 
100 time steps were integrated from the initial redshift 2; = 25 to Zmed- 

We then select random lines-of-sight through the simulation box along which to 
interpolate one-dimensional density and velocity field using the Daubechies-4 scaling 
functions with j=7. Using one-dimensional density field, we assign temperature to each 
pixel using the polytropic equation of state T = Top", where Tq and a depend on the 
spectral shape of the UV background and on the history of reionization. We adopt typical 
values To ~ lO^K and a = 0.6. The neutral hydrogen fraction in each pixel is computed 
by adopting the cosmic abundance of hydrogen and assuming photoionization equilibrium. 
The optical depth r at a given pixel is then obtained by integrating in real space by 
including the effect of peculiar velocity field and convolving with Voigt thermal broadening 
profile. To match the parameters of observed spectra, r was computed onto a 2^ = 512 grid 
and smoothed by a gaussian window to match with the spectral resolution of observations. 
The absorption transmitted spectra F = exp(— r) is normalized such that the mean 
transmission in the spectra matches with observations. We take 29 different values of the 
mean transmission, at different redshifts. 



We calculated the effective mean transmission e~'^^^'^j [eq.(49)] for the pseudo- hydro 
simulation samples. As an example, we show a result of e~^('*')j vs. j for a simulated sample 
in Fig. 5. It is similar to Fig. 3, i.e. the effective mean transmission on small scales is 
lower than the mean transmission. Fig. 6 gives 29 e^'^^^'^j on scale k ~ 0.02 h Mpc~^, and 
comparing them with the 29 e^^ used for the simulation. It shows that all the effective 
mean transmission on this small scale is lower than the corresponding mean transmission 
used for simulation. This result is the same as lognormal model. 

Therefore, the self-normalized power spectrum estimator eq.(41) is sensitive to both the 
mean transmission (on large scales) and the anti-correlation between the fluctuation and 
background on small scales. The power spectrum given by traditional normalization is only 
sensitive to the mean transmission, but not the fluctuation-background anti-correlation on 
small scales. Therefore, the estimator eq.(41) would be more useful for the discrimination 
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between models of Lya forests. We will demonstrate this point in next subsection. 



5.3. The power spectrum of a real sample 

To demonstrate the estimator eq.(41), we now analyze a real data set: the Lya forest 
spectrum of ql700+642 given by Dobrzycki and Bechtold (1996). The data include the flux 
and continuum in the wavelength range from 3731.04 to 4611.26 A with resolution ~ 0.25 
A. We use the 2^^ pixels from 3816.81 to 4080.38 A. We remove regions in the spectra where 
the signal to noise ratio is less than 2.0. 

Using the estimator in eq.(41), we calculated the normalized DWT power spectrum 
Pj. From the data of flux and continuum, one can also calculate the unnormalized power 
spectrum Pj by eq.(46). Figure 7 shows the power spectra Pj and P*. It also shows the 
traditionally normalized power spectrum, i.e. Pj/e~^. 

Similar to Fig. 2, Fig. 7 shows that the shapes of the two power spectra Pj and Pj 
of ql700+642 are completely the same on large scales (j < 6, or > 1.1 h^^ Mpc). In this 
scale range, the power spectrum Pj is perfectly the same as the traditionally normalized 
power spectrum (see the coincidence of the solid lines and dot-dashed line). It means 
that the estimator eq.(41) gives correct normalization on large scales. The ratio of Pj/ Pj 
on large scales (j < 6, or > 1.1 h~^ Mpc) is ~ 1.4, and therefore the mean transmission 
is e~'^ ~ I/a/1.4 = 0.845. This value is exactly the same as the e""^ given by direct 
measurement. 

On small scales j > 6 (1.1 h^^ Mpc), the anti-correlation between the background 
and fluctuations becomes signiflcant. Fig. 8 shows that the WFC- SFC anti-correlation is 
stronger on scales j = 8 and 9 (100 - 300 h^^ kpc). In this scale range, the self-normalized 
Pj is generally larger than that given by traditional normalization, i.e. with a constant 
mean transmission. At j = 9 (~ 140 h^^ kpc), the traditionally normalized power spectrum 
is lower than the self-normalized Pj by a factor of 5.6. 

What does this difference indicate? We know that the WFC-SFC anti-correlation is 
due to the intermittency of the fleld. In terms of Fourier decomposition, an intermittent 
fleld means that the phases of the Fourier modes is highly correlated, but the one point 
distribution of the Fourier amplitude is still gaussian due to the central limit theorem. 
Therefore, the WFC-SFC anti-correlation can be eliminated by a phase-randomization, 
which is produced by taking the inverse transform of the Fourier coefficients of the 
original data after randomizing their phases uniformly over [0, 27r] without changing their 
amplitudes. The phase randomized sample gets rid of the intermittent behavior possessed 
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by the field, but the unnormahzed power spectrum P* and the traditional normalized 
power spectrum of the phase randomized sample are exactly the same as the original 
one (Jamkhedkar, Zhan & Fang 2000; Zhan, Jamkhedkar & Fang 2001). In Fig. 7, the 
long dashed line is the unnormahzed power spectrum of the original data and its phase 
randomized counterpart, and the dot dashed line is the traditionally normalized power 
spectrum of the original data and its phase randomized counterpart. 

Therefore, neither unnormahzed nor traditionally normalized power spectrum can 
distinguish between the highly intermittent field and its phase-randomized counterpart. 

On the other hand, the estimator eq.(41) can detect the difference between the 
two fields. We calculate the power spectrum of the phase randomized sample by the 
self- normalized estimator eq.(41). The result is shown in Fig. 7. It shows that the self- 
normalized power spectrum of the original data is very different from its phase-randomized 
counterpart. 

Therefore, one can conclude that the self-normalized power spectrum estimator is 
sensitive to both the clustering behaviors of the field on large scales (mean transmission) 
and on small scales (intermittency, or fluctuation-background anti-correlation). But the 
traditional normalized power spectrum is insensitive to the phase correlation of the Fourier 
modes. 



6. Conclusion 

The power spectrum of Lya forests is a direct indicator of the matter distribution 
at high redshift. This paper addresses the issue of how continuum fitting and the mean 
transmission affect the estimated power spectrum of QSO's Lya forests. We propose a 
straightforward method for calculating the power spectrum of observed Lya forests. This 
method is based on the DWT decomposition of the transmission flux. It gives a consistent 
calculation for the decomposition of flux and the normalization of power spectrum. 

With numerical simulation samples, we showed that the power spectrum obtained by 
this estimator is independent of the continuum. The non-linear power spectrum of the 
transmission can be reliably recovered from the observed flux regardless of the continuum, 
i.e. the algorithm can automatically take care of the normalization by the continuum 
without a continuum fitting. 

With numerical simulation samples, we also show that the power spectrum estimator 
can automatically consider the normalization of the mean transmission, i.e. the algorithm 
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doesn't need a pre-calculated mean transmission to do the normalization. 

For a gaussian field, the power spectrum given by the proposed estimator principally 
is the same as the power spectrum given by traditional normalization. In this case, an 
advantage of the proposed estimator is that it is free from fitting parameters. 

On scales with significant non-linear clustering, like intermittency or phase correlation, 
the self-normalized power spectrum is essentially different from the power spectrum 
normalized by traditional method. The latter is not sensitive to the phase correlation, while 
the former is. Therefore, as an estimator of power spectrum of non-linear field traced by 
Lya forests, the self-normalization algorithm is useful of the discrimination among models. 

We thank Drs. L.L. Feng and W.L.Lee for their help. PJ would also like to thank 
Jennifer Scott for useful discussions. 



A. Photon counts and flux 

The observed photon counts, N{X), at a pixel corresponding to wavelength A is given 
by 

Nob{X) = C(A)iV,(A)e-^W + E{X) (Al) 

where r(A) is the optical depth, and Nc is the continuum. C(A) describes the A-dependence 
of CCD's efficiency. E[X) is noise, whose mean and covariance are 

{E{X)) = E{X), {E{X)E{X')) = 45f,„ (A2) 

where 6^^ is Kronecker delta function. Eq.(A2) means that the random variable E{X) is 
independent for each wavelength A. 

The observed photon counts are reduced as 

^,(,) . MA1_M. (A3) 

We have then 

A(A) = A,(A)e-^W - ^^[^(A) - i?(A)]. (A4) 

If we define a new variable for error as 

N,{X) = -^^[E{X) ~ E{X)], (A5) 
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we have 



where 



iV(A)=iV,(A)e--W-iV,(A), 



(iVe(A)) = 0, (iVe(A)iVe(A')) = cr"{X)6^^,, 

The fluctuation of the transmission is defined by 



5{\) 



rW _ r(A) 



-r(A) 



Thus, eq.(A6) yields eq.(21). 



(A6) 

(A7) 

(A8) 



B. Poisson sampling and modified Poisson sampling 

Consider the reduced photon count N{x) as a sampling of random field 

N{x) =N{x)[l + S{x)], 



(Bl) 



where N{x) = Net '^(''*). For the Poisson sampling, the characteristic function of N{x) is 

2^^r^N{x)u(x)dx^ = exp I /rfxiV(x)[e*"(") - 1]| . (B2) 



Thus, the correlation functions of N{x) are given by 



{N{x,)...N{xn))p = ^ 



S^'Z 



6u{xi)...6u{xn) 



u=0 



where {■■■)p is the average for the Poisson sampling. We have then 

{N)p = F{x), 
and 



{N{x)N{x'))p = N{x)N{x') + 6^{x - x')N{x). 



This equation yields 



^^^^'^^^^'■>^-^^{Ww)-'"^'-'^m 



This gives eq.(24). 



(B3) 

(B4) 
(B5) 

(B6) 
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For a weighted Poisson sampling, the data at x are given as a Poisson samphng of 
N{x), but with a weight g{x). In this case, the characteristic function eq.(B2) becomes 

2^^ijpnx)u(x)dx^ = exp I /"rfa;iV(a;)[e*3(^)"(^) - 1]| . (B7) 

We have then 

{N)p = g{x)N{x), (B8) 

and 

{N{x)N{x'))p = g{x)g{x')N{x)N{x') + 6^{x - x')g\x)N{x). (B9) 

Eq.(B2) yields 

(5(x)5(x')) = -1 + l^W^^T^xl^) - ^"(- - -')t^- (BIO) 

^ ^ ' ^ " \g{x)N{x)g{x)N{x') ^ ' N{x) ^ ' 
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Fig. 1. — A simulation sample of the QSO Lya forests. The transmission e ^ is shown in 
top panel. The lower four panels are the Poisson sampling of flux F^.e^'^ plus noise with 1.) 
F, = 2000; 2.) F^ = 2000(A/Ao)^ Aq = 41331, 3.) F, = 2000[1 + 0.1 sin[27r(A-A„i„)/(A„ax- 
Xmin)]; and 4.) Fc = 2000[1 + 0.5 sin[27r(A - Xmin)/{>^max - Xmin)], respectively. 

Fig. 2. — Upper panel: The DWT power spectra Pj of the simulated Lya forest samples 
with continuum as 1.) F^ = 2000 (solid line); 2.) F^ = 2000(A/Ao)^ (dotted hue), where 
Aq = 4133 A; and 3.) F^. = 2000[1 + 0.1 sin[27r(A — \min)/iXmax ~ ^min)] (short-dotted line), 
respectively. Lower panel: The DWT power spectra Pj of the simulated Lya forest samples 
with continua given as Fc = 2000[1 + (isin[27r(A — \min)/{Xmax — Amm)] and d = 0.1 (sohd), 
0.2 (dotted), 0.3(short dashed), 0.4 (long dashed) and 0.5 (dot-dashed line), respectively. 
The short-dash long-dashed line is the power spectrum of the transmission, P*. The 1-sigma 
error bars are calculated from 100 samples. The scale j corresponds to 2^^~^ x 2 h~^ kpc. 



Fig. 3. — The effective mean transmission e~'^^'^^j vs. j (circle) [eq.(49)] for three lognormal 
simulation samples. The solid line is the mean transmission used for the simulation. 

Fig. 4. — The histogram showing the number distribution of the SFCs. Their values can be 
read from the right hand y axis. The squares show the average value of the — WFCs — in 
each SFC bin. Their values can be read from the left hand y axis. The scale j corresponds 
to 2^^-^ X 2 h-^ kpc. 



Fig. 5. — The effective mean transmission e '^(^^j vs. j (circle) [eq.(49)] for a pseudo-hydro 
simulation samples. The solid line is the mean transmission used for the simulation. 



Fig. 6. — The effective mean transmission e '^^^^j = JPj/Pj vs. the mean transmission e "^ 
used for the 29 pseudo-hydro simulations. The solid line is for e""^. 

Fig. 7. — The transmission power spectrum Pj (solid line) for the real sample ql 700-1-642, 
calculated by estimator (41), and the power spectrum Pj (long dashed line) calculated 
by eq.(46). The dot-dashed line represents the power spectrum normalized by traditional 
method. The dotted line is the power spectrum by eq.(46) for the phase randomized field. 
The error bar is given by the maximum and minimum of 20 realizations. The scale j 
corresponds to 2^^^^ x 34 h^^ kpc. 

Fig. 8. — The WFC-SFC anti-correlation for real sample of the ql 700-1-642 Lya forests. The 
histogram showing the number distribution of the SFCs. Their values can be read from the 
right hand y axis. The squares show the average value of the WFCs in each SFC bin. Their 
values can be read from the left hand y axis. The scale j corresponds to 2^^~^ x 34 h~^ kpc. 
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